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This  paper  proposes  a  method  for  short  term  security-constrained  unit  commitment  (SCUC)  for  hydro 
and  thermal  generation  units.  The  SCUC  problem  is  modeled  as  a  multi-objective  problem  to 
concurrently  minimize  the  ISO's  cost  as  well  as  minimizing  the  emissions  caused  by  thermal  units. 
The  non-linearity  of  valve  loading  effects  is  linearized  in  the  presented  problem.  In  order  to  model  the 
SCUC  problem  more  realistically,  this  paper  considers  the  dynamic  ramp  rate  of  thermal  units  instead  of 
the  fixed  rate.  Moreover,  multi-performance  curves  pertaining  to  hydro  units  are  developed  and  the 
proposed  SCUC  problem  includes  the  prohibited  operating  zones  (POZs).  Besides,  the  model  of  SCUC  is 
transformed  into  mixed  integer  linear  programming  (MILP)  instead  of  using  mixed  integer  non-linear 
programming  (MINLP)  which  has  the  capability  to  be  solved  efficiently  using  optimization  software  even 
for  real  size  power  systems.  Pareto  optimal  solutions  are  generated  by  employing  lexicographic 
optimization  as  well  as  hybrid  augmented-weighted  e-constraint  technique.  Furthermore,  a  Fuzzy 
decision  maker  is  utilized  in  this  paper  to  determine  the  most  preferred  solution  among  Pareto  optimal 
solutions  derived  through  solving  the  proposed  multi-objective  SCUC  problem.  Eventually,  the  proposed 
model  is  implemented  on  modified  IEEE  118-bus  system  comprising  54  thermal  units  and  8  hydro  units. 
The  simulation  results  reveal  that  the  solutions  obtained  from  the  proposed  technique  in  comparison 
with  other  methods  established  recently  are  superior  in  terms  of  total  cost  and  emission  output. 
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1.  Introduction 

In  vertically  integrated  power  systems,  the  goal  of  implementing 
unit  commitment  (UC)  is  optimizing  energy  sources  in  a  way  that  the 
load  demand  is  supplied  with  the  least  cost.  Besides,  if  security  criteria 
are  included  in  the  UC,  then  it  is  known  as  security-constrained  unit 
commitment  (SCUC).  If  the  UC  is  performed  by  each  of  the  generation 
companies  (GENCOs)  to  optimize  energy  sources  in  a  way  that  the 
GENCO’s  profit  becomes  maximized,  it  is  called  Price-Based  Unit 
Commitment  (PBUC)  [1].  A  MILP  framework  for  the  UC  problem 
relating  to  thermal  units  is  proposed  in  Ref.  [2]  in  which  fewer  number 
of  binary  variables  and  constraints  are  needed  compared  to  the 
models  reviewed  before  leading  to  remarkable  less  computational 
burden.  On  the  basis  of  the  SCUC  model  presented  in  Ref.  [3],  an 
efficient  AC  corrective/preventive  contingency  dispatch  comprising 
24-h  time  period  is  proposed.  A  SCUC  model  considering  energy  and 
ancillary  sen/ices  auction  is  proposed  in  Ref.  [4]  which  is  applicable  to 
reserve  requirements  optimization  in  electricity  market  framework  by 
ISO.  Ref.  [5]  has  suggested  a  model  having  the  capability  to  solve  the 
coordinated  generation  and  transmission  maintenance  scheduling 
with  SCUC  in  which  the  scheduling  period  is  weeks  to  months  and 
decomposition  of  the  optimization  problem  is  performed  using 
Lagrangian  relaxation  method.  In  Ref.  [6],  the  intermittent  and  vola¬ 
tile  wind  power  generation  is  included  in  the  SCUC  problem 
while  Benders  decomposition  is  employed  to  solve  the  problem.  The 
problem  of  determining  optimal  response  of  a  thermal  unit  along  with 
a  hydro  GENCO  to  a  spot  electricity  market  is  investigated  in  Refs.  [7,8] 
in  which  the  profit  of  the  unit  through  participating  in  spot  market  to 
sell  both  energy  and  spinning  reserve  is  set  as  the  objective  function  to 
be  maximized.  The  problem  of  short-  hydrothermal  scheduling  is 
addressed  in  Ref.  [9]  utilizing  Modified  Differential  Evolution  (MDE) 
algorithm.  Ref.  [10]  has  proposed  mid-term  risk-constrained  hydro- 
thermal  scheduling  problem  in  a  stochastic  framework  for  a  GENCO 
that  its  objective  is  assigned  to  maximize  payoffs  as  well  as  minimiz¬ 
ing  financial  risks  during  its  mid-term  scheduling  of  thermal,  cascaded 
hydro  in  addition  to  pumped-storage  units.  The  act  of  determining  the 
PBUC  status  for  power  generation  companies  is  performed  in  Ref.  [10] 
through  solving  self-scheduling  problem  prior  to  giving  their  bids  to  a 
day-ahead  market.  The  proposed  model  is  formulated  as  a  determi¬ 
nistic  optimization  problem  to  maximize  the  profit  exploiting  0/1 
MILP  technique.  Ref.  [11]  has  employed  Improved  Particle  Swarm 
Optimization  (IPSO)  algorithm  to  solve  the  short-term  hydrothermal 
scheduling  problem  in  order  to  obtain  the  optimal  power  generation. 
The  GENCO's  arbitrage  problem  is  taken  into  consideration  in  Ref.  [12] 
utilizing  stochastic  PBUC  while  modeling  the  associated  risks.  The 
requirements  of  short  term  hydro-thermal  scheduling  (SHTS)  problem 
are  not  met  anymore  after  passage  of  the  clean  air  act  amendments  in 
1990  [13]  and  the  emission  issues  must  be  unavoidably  taken  into 
account.  Ref.  [14]  was  one  of  the  foremost  work  trying  to  solve  the  UC 
problem  while  emissions  minimization  had  been  taken  into 


consideration  as  another  objective  function  and  the  proposed  problem 
was  solved  employing  weighting  method  to  make  the  two-objective 
problem  a  single-objective  optimization  problem.  Unsupported  effi¬ 
cient  solutions  cannot  be  generated  by  weighting  method  in  multi¬ 
objective  integer  problems  as  well  as  MIP  ones.  The  results  derived 
from  weighting  method  are  highly  dependent  on  the  scaling  of 
objective  functions.  Thus,  objective  functions  should  be  scaled  to  a 
common  scale  prior  to  the  weighted  sum  formation  [15].  Multi¬ 
objective  evolutionary  algorithm  is  applied  to  the  environmental/ 
economic  dispatch  (EED)  problem  in  Ref.  [IS]  to  simultaneously 
minimize  the  fuel  cost  as  well  as  minimizing  the  emissions.  Unlike 
weighting  method  needing  large  number  of  optimization  problems  to 
be  solved,  the  suggested  approach  in  Ref.  [16]  requires  only  a  single 
run  to  trace  the  Pareto  optimal  front.  Ref.  [17]  solves  the  stochastic 
multi-objective  generation  scheduling  problem  employing  weight  age 
pattern  search  methods  while  operating  cost,  NOx  and  the  risk  caused 
by  the  variance  of  active  as  well  as  reactive  power  mismatch  are  set  to 
be  minimized  concurrently  and  the  proposed  model  is  implemented 
on  IEEE  11 -bus  system.  The  problem  of  minimizing  cost,  S02,  C02,  and 
NOx  is  addressed  in  Ref.  [18]  using  Non-Inferior  Surface  estimation 
method.  The  electric  power  dispatch  is  done  in  Ref.  [19]  employing  a 
fuzzified  multi-objective  particle  swarm  optimization  (MOPSO)  algo¬ 
rithm  while  economic  and  environmental  issues  are  all  taken  into 
consideration  and  the  results  derived  from  this  technique  is  compared 
with  those  obtained  from  weighted  sum  method  and  evolutionary 
multi-objective  optimization  algorithms.  Ref.  [20]  presents  an  inter¬ 
active  Fuzzy  satisfying  approach  on  the  basis  of  evolutionary  program¬ 
ming  (EP)  technique  while  the  decision  maker  (DM)  is  supposed  to 
have  Fuzzy  targets  for  each  objective  function;  afterwards  the  corre¬ 
sponding  Pareto  optimal  solution  for  the  decision  maker  targets  is 
generated  by  implementing  the  Fuzzy  satisfying  method  based  on  the 
EP  technique.  A  stochastic  multi-objective  framework  is  presented  in 
Ref.  [21]  utilizing  Fuzzy  decision  making  wherein  the  tradeoff  relation 
between  conflicting  objective  functions  is  simulated  by  weighting 
method,  e-constraint  is  another  multi-objective  technique  utilized  in 
Ref.  [22]  to  produce  non-dominated  solutions  for  thermal  power 
dispatch  with  operating  cost  and  emissions  as  objective  functions  that 
are  assumed  to  be  minimized;  After  that,  the  most  preferred  solution 
is  determined  for  the  proposed  multi-objective  problem  by  surrogate 
worth  trade-off  approach.  In  Refs.  [23,24],  the  economic/emission  load 
dispatch  problem  is  formulated  in  the  multi-objective  framework  and 
solved  employing  Hopfield  neural  networks  and  non-dominated 
sorting  genetic  algorithm  (NSGA),  respectively.  Ref.  [25]  suggests 
MOPSO  algorithm  to  solve  the  stochastic  economic/emission  load 
dispatch  problem.  The  NSGA  is  utilized  in  Ref.  [26]  to  generate  the 
Pareto  optimal  front  and  the  most  compromise  solution  is  selected  by 
multi-attribute  decision  making  method.  The  e-constraint  is  applied  to 
the  multi-objective  problem  with  objective  functions  of  GENCO's  profit 
to  be  minimized  and  emissions  to  be  minimized  in  a  day-ahead  joint 
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Nomenclature 

Indices 

I  index  for  thermal  units 

j  index  for  hydro  units 

f  index  for  time  interval  (h) 

Constants 

;/  conversion  factor  equal  to  3.6  x  10~3  (hm3/s/m3/h) 

0  time  intervals  of  the  considered  planning  horizon 

0(j,  t)  minimum  water  discharge  of  hydro  unit  j  at  time 
t  (m3/s) 

0(j,  t )  maximum  water  discharge  of  hydro  unit  j  at  time 

t  (m3/s) 

t,j  time  delay  of  water  transport  between  hydro  unit 

i  reservoir  and  hydro  unit  j  reservoir  (h) 

A,  shut-down  cost  of  thermal  unit  i  ($) 

Aj  start-up  cost  of  hydro  unit  j  ($) 

bn(i)  slope  of  block  n  of  fuel  cost  curve  pertaining  to 

thermal  unit  i  ($/MW  h) 

bn(j)  slope  of  the  volume  block  n  of  the  reservoir  pertaining 
to  hydro  unit./'  (m3/s/hm3)  (1  hm3  =  106  m3) 
bkn(j)  slope  of  the  block  n  pertaining  to  the  performance 
curve  k  of  hydro  unit  j  (MW/m3/s) 
ben(i)  slope  of  segment  n  in  emission  curve  of  thermal  unit 
i  (lbs/MW  h) 

Cj  valve  loading  coefficient 

emin(0  minimum  emission  generation  of  thermal  unit  i  (lbs) 
E(p“ _1(i))  emission  generation  of  n-lth  upper  limit  in  emis¬ 
sion  curve  of  thermal  unit  i  (lbs) 
ft  valve  loading  coefficient 

F(p"_ ,(!))  cost  of  generation  of  n-lth  upper  limit  in  fuel  cost 
curve  pertaining  to  the  thermal  unit  i($/h) 

F(J,t)  forecasted  natural  water  inflow  of  the  reservoir  asso¬ 
ciated  to  the  hydro  unit  j  at  hour  t  (hm3/h) 

L  number  of  variable  head 

M  number  of  prohibited  operating  zones 

NL  number  of  piecewise  linearization  blocks  pertaining  to 

the  variable  cost  function 

minimum  output  power  of  thermal  unit  i  (MW) 
Pmax{i )  maximum  output  power  of  thermal  unit  i  (MW) 

Pn0)  minimum  output  power  of  hydro  unit  j  for  perfor¬ 

mance  curve  n  (MW) 
p(j)  capacity  of  hydro  unit  j  (MW) 

pd(i)  lower  limit  of  nth  prohibited  operating  zone  pertain¬ 
ing  to  thermal  unit  i  (MW) 

p“_  ,(i)  upper  limit  of  n-lth  prohibited  operating  zone  per¬ 
taining  to  thermal  unit  i  (MW) 

Q_(j)  minimum  water  discharge  of  hydro  unit  j  provided 
that  it  is  on  (m3/s) 

Q_n(j)  maximum  water  discharge  of  block  n  of  hydro  unit 
j  (m3/s) 

RDLn(i )  ramp-down  limit  for  block  n  (MW) 

RULn(i)  ramp-up  limit  for  block  n  (MW) 

SUE(i)  emission  generated  by  thermal  unit  i  when  started- 
up  (lbs) 

SDE(i)  emission  generated  by  thermal  unit  i  when  shut¬ 
down  (lbs) 

SD(i)  shut-down  ramp  rate  limit  of  thermal  unit  i  (MW  h) 
SU(i)  start-up  ramp  rate  limit  of  thermal  unit  i  (MW  h) 
v0(j )  minimum  volume  water  of  the  reservoir  pertaining  to 

hydro  unit  j  (hm3) 


v°(j)  volume  of  the  reservoir  at  the  beginning  of  the 
planning  horizon  (hm3) 

v8(j)  volume  of  the  reservoir  at  the  end  of  the  planning 
horizon  (hm3) 

v„0)  maximum  volume  of  the  reservoir  j  pertaining  to  the 
nth  variable  head  (hm3) 

Variables 


pn(i,  t )  binary  variable  that  is  equal  to  1  if  block  n  of  fuel  cost 

curve  of  thermal  unit  i  at  time  t  is  selected 
/i„0,  t)  binary  variable  that  is  equal  to  1  provided  that  the 
variable  head  n+1  of  hydro  unit  j  at  time  t  is  selected 
X„ (i,  t )  binary  variable  that  is  equal  to  1  provided  that  the 

output  power  of  thermal  unit  i  at  time  t  has  over¬ 
stepped  block  n  of  valve  loading  effects  curve 
S„(i,  t)  generation  of  block  n  of  fuel  cost  curve  of  thermal  unit 
i  at  time  t  (MW) 

<Z>  payoff  table 

i//  n(i,  t )  generation  of  block  n  pertaining  to  the  thermal  unit  i 
at  time  t  of  valve  loading  effects  curve  (MW) 
t)  volume  of  block  n  for  the  reservoir  of  hydro  unit  j  at 
time  t  (hm3) 

£2  feasible  region 

/  total  membership  function  of  the  rth  Pareto  optimal 

solution 

degree  of  optimality  for  the  nth  objective  function  in 
the  rth  Pareto  optimal  solution 
B(i,t)  start-up  cost  of  thermal  unit  i  at  time  t  ($) 

C(i,t)  cost  of  valve  loading  effects  pertaining  to  the  thermal 
unit  i  at  time  t  ($) 

F(i,t)  fuel  cost  of  thermal  unit  i  at  time  t  ($) 

Fcost  main  objective  function 

^mission  seconcj  objective  function 

hn (j,t)  binary  variable  that  is  equal  to  1  provided  that  the 
water  discharge  of  hydro  unit  j  at  time  t  has  over¬ 
stepped  block  n 

I(i,t)  binary  variable  that  is  equal  to  1  provided  that  the 
thermal  unit  i  is  on-line  at  time  t 
I(j,t)  binary  variable  that  is  equal  to  1  provided  that  the 
hydro  unit  j  is  on-line  at  time  t 

Id(i,t)  binary  variable  that  is  equal  to  1  if  the  thermal  unit  i  at 
time  t  provides  non-spinning  reserve  when  unit  is  off 
P(i,t)  output  power  of  thermal  unit  i  at  time  t  (MW) 

p(  i,  t)  maximum  output  power  of  thermal  unit  i  at  time  t  (MW) 

p(j,t )  output  power  of  hydro  unit  j  at  time  t  (MW) 

Q(j,t)  water  discharge  of  hydro  unit  j  at  time  t  (m3/s) 
qn(j !,t)  water  discharge  of  block  n  pertaining  to  hydro  unit  j  at 

time  t  (m3/s) 

RDL(p(i,t ))  ramping-down  limit  of  thermal  unit  i  at  time  t  (MW) 
RUL(p(i,t ))  ramping-up  limit  of  thermal  unit  i  at  time  t  (MW) 
s(j,t )  spillage  of  the  reservoir  associated  to  hydro  unit  j  at 

hour  t  (m3/s) 

v(j,t)  amount  of  water  pertaining  to  hydro  unit  j  reservoir  at 

hour  t  (hm3) 

y(i,t)  binary  variable  that  is  equal  to  1  provided  that 

thermal  unit  i  is  started-up  at  the  beginning  of  time  t 
y(J,t)  binary  variable  that  is  equal  to  1  provided  that  hydro 
unit  j  is  started-up  at  the  beginning  of  time  t 
z(i,t )  binary  variable  that  is  equal  to  1  provided  that 

thermal  unit  i  is  shut-down  at  the  beginning  of  time  t 
z(j,t )  binary  variable  that  is  equal  to  1  provided  that  hydro 

unit  j  is  shut-down  at  the  beginning  of  time  t 
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Sets 

I  set  of  thermal  generating  units 

J  set  of  hydro  generating  units 

N  set  of  indices  of  blocks  of  the  piecewise  linearization 

of  hydro  units’  performance  curve 


NE  set  of  indices  of  blocks  of  the  piecewise  linearization 
of  thermal  units'  emission  curve 
T  set  of  indices  of  periods  of  the  market  time  horizon 

£2j  set  of  upstream  reservoirs  of  hydro  unit  j 


reserve  and  energy  market.  The  problem  of  short-term  multi-objective 
hydrothermal  scheduling  taking  into  consideration  the  cost  and 
emissions  as  two  objective  functions  is  solved  using  in  Refs.  [29-31] 
using  differential  evolution,  quantum-behaved  particle  swarm 
optimization,  and  quantum-behaved  particle  swarm  optimization 
techniques,  respectively  while  ignoring  security  issues  as  well  as 
implementing  the  model  on  a  case  study  with  small  size.  Further 
information  on  optimization  techniques  presented  to  solve  the  SHTS 
problem  can  be  found  in  Refs.  [32,33], 

This  paper  solves  the  multi-objective  problem  with  total  cost  and 
emission  as  objective  functions  that  are  intended  to  be  minimized 
employing  lexicographic  optimization  and  augmented  e-constraint 
technique.  This  approach  can  be  well  implemented  by  the  ISO  for 
scheduling  wherein  more  practical  constraints  of  thermal  and  hydro 
units  are  considered  in  comparison  with  previous  works  in  this  area. 
So  far,  there  is  not  any  work  devoted  to  solve  the  SCUC  problem 
minimizing  total  cost  as  well  as  minimizing  emissions  for  a  large- 
scale  case  study  (IEEE  118-bus  system)  taking  into  consideration  high 
number  of  practical  limitations,  such  as  Prohibited  Operation  Zones 
(POZs),  variable  start-up  cost  function,  minimum  up-time  (MUT)  and 
minimum  down-time  (MDT),  valve  loading  effects,  dynamic  ramp-up 
limit  (RUL)  and  ramp-down  limit  (RDL)  for  thermal  units  as  well  as 
multi-performance  curves  for  hydro  units.  The  main  contributions  of 
this  paper  are  listed  below: 

•  Proposing  a  multi-objective  framework  for  the  SCUC  problem 
taking  cost  and  emissions  as  the  two  objective  functions. 

•  Employing  lexicographic  optimization  along  with  augmented 
e-constraint  method  to  simultaneously  minimize  the  total  cost 
and  emissions  of  thermal  generating  units. 

•  Presenting  valve  loading  effect  with  linear  formulation  as  well 
as  replacing  fixed  ramp  rate  limit  with  dynamic  one  for 
thermal  units. 

•  Using  flexible  approach  for  multi  POZs  of  thermal  generating 
units  as  well  as  multi-performance  curves  for  hydro  units. 

•  Linearizing  the  MINLP  problem  and  converting  it  to  MILP 
problem  and  using  an  effective  analytical  technique  to  solve 
the  proposed  problem. 

The  remainder  of  this  paper  is  categorized  as  below: 

The  proposed  lexicographic  optimization  and  augmented  E-con¬ 
straint  method  are  included  in  Section  2.  The  thermal  and  hydro 
model  is  presented  in  Section  3  while  Section  4  comprises  the  case 
studies  and  represents  the  results  with  detailed  discussion.  Finally, 
some  relevant  conclusions  of  this  paper  are  drawn  in  Section  5. 


optimality  i.e.  the  solution  that  cannot  get  better  except  making 
its  performance  worse  in  at  least  one  of  the  other  objective 
functions.  An  effective  approach  presented  to  deal  with  multi¬ 
objective  problems  is  s-constraint  technique  having  one  main 
objective  function  selected  among  all  objective  functions.  The 
E-constraint  method  is  prior  to  traditional  weighting  methods 
from  several  aspects.  The  deficiency  by  the  weighting  methods 
is  forming  a  single  objective  function  out  of  multi-objective 
problem  by  uniting  the  objective  functions  using  weighted  sum. 
The  advantages  of  utilizing  E-constraint  can  be  listed  as  below 
[15,27,34]: 

•  The  weighting  method  results  in  merely  efficient  extreme 
solutions  for  linear  problems  whereas  non-extreme  solutions 
can  be  generated  through  applying  s-constraint  approach. 

•  Unlike  the  weighting  method,  in  multi-objective  integer  pro¬ 
gramming  problems  as  well  as  MIP  ones,  unsupported  efficient 
solutions  can  be  generated  by  £-constraint. 

•  The  weighting  method  highly  depends  on  the  scaling  of 
objective  functions  while  this  issue  is  not  important  in  the 
E-constraint. 

•  The  controlled  number  of  efficient  solutions  can  be  generated 
in  £-constraint  technique  by  effectively  tuning  the  number  of 
grid  points  in  the  range  of  each  objective  function. 

In  spite  of  the  advantages  of  E-constraint  approach  have,  there  are 
two  key  points  that  should  be  noted.  First,  there  is  no  guarantee  that 
the  range  of  the  objective  functions  over  the  efficient  set  is  optimal 
while  lexicographic  optimization  method  is  employed  to  overcome 
this  deficit.  The  second  problem  with  e-constraint  is  the  possibility  of 
being  dominated  or  inefficient  of  the  produced  Pareto  optimal 
solutions.  One  likely  way  to  eliminate  this  shortfall  is  to  use  augmen¬ 
ted  E-constraint  method.  Further  details  on  employing  lexicographic 
optimization  and  augmented  £-constraint  can  be  found  in  the  previous 
paper  of  the  authors  in  this  area  [34]  while  this  procedure  is  not 
repeated  in  this  paper.  The  importance  of  the  objective  functions  is  not 
taken  into  consideration  by  the  augmented  E-constraint  method 
during  Pareto  solutions  generation.  On  the  contrary,  the  augmented- 
weighted  E-constraint  method  clearly  considers  the  relative  impor¬ 
tance  of  the  objective  functions  in  the  Pareto  solutions  generation. 
In  this  paper,  the  MMP  model  of  SCUC  comprising  two  objective 
functions  i.e.  Fh  F2  is  solved  using  lexicographic  optimization  and 
hybrid  augmented-weighted  £-constraint  method.  The  detailed 
description  of  the  MMP  solution  methodology  is  given  in  the  next 
section. 


2.  Multi-objective  mathematical  programming 

The  concept  of  multi-objective  mathematical  programming 
(MMP)  involves  with  more  than  one  objective  function  and 
there  is  no  longer  a  single  optimal  solution  that  concurrently 
optimizes  all  objective  functions.  In  such  situations,  the  most 
desired  solution  is  intended  to  be  derived  by  the  DM.  Also,  the 
concept  of  optimality  is  substituted  with  efficiency  or  Pareto 


2.3.  Solution  methodology  of  MMP 


The  MMP  formulation  on  the  basis  of  the  augmented  s-con- 
straint  method  is  presented  as  (1)  while  Refs.  [27,34]  include 
further  details: 


Min  /  Max 


Fj(x)+din  r, 


s.t.  Fi(x)-dirjSj  =  e,  V  s,-  e  R+ ,  i  =  2,3,  ...,p. 


(1) 
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where,  the  number  of  objective  functions  is  denoted  by  p.  dir, 
expresses  the  direction  of  the  ith  objective  function  that  in  the 
case  of  minimizing  the  objective  function,  dir,-  is  - 1  while  + 1 
shows  that  the  objective  function  is  supposed  to  be  maximized. 
The  efficient  solutions  of  the  problem  are  derived  by  parametrical 
iterative  variations  in  et.  The  introduced  slack  or  surplus  variables 
for  the  MMP  problem's  constraints  is  denoted  by  s,.  Moreover,  the 
second  term  of  the  MMP  problem  comprises  ri(s,-/r,-)  to  prevent 
any  scaling  problem.  Because  of  the  augmentation  of  the  objective 
function  F3  by  the  second  term,  formulation  (1)  is  called  augmen¬ 
ted  e-constraint  approach.  The  generation  of  only  efficient  solu¬ 
tions  by  the  augmented  e-constraint  can  be  expressed  [34]  while  it 
is  proven  in  Ref.  [15],  The  ranges  of  objective  functions  should  be 
determined  to  implement  the  e-constraint  approach.  These  ranges 
can  be  computed  from  payoff  table  as  the  most  usual  way.  The 
ordinary  approach  to  form  the  payoff  table  does  not  ensure  that 
the  solutions  derived  from  individual  optimization  of  the  objective 
functions  are  efficient  solutions  or  Pareto  optima.  Hence,  lexico¬ 
graphic  optimization  is  exploited  in  this  paper  to  eliminate  the 
defect  by  inefficient  solutions  once  the  payoff  table  is  formed. 
In  general,  the  procedure  of  lexicographic  optimization  of  several 
objective  functions  is  to  take  and  optimize  the  first  objective 
function  while  between  possible  alternative  optima  optimize  for 
the  second  objective  function  and  so  on.  Ref.  [34]  includes  the 
detailed  procedure  of  calculating  payoff  table  for  the  MMP 
problem.  The  Payoff  table  is  a  p  by  p  table.  The  obtained  values 
of  the  objective  function  F,-  are  included  in  the  ith  column  of 
the  payoff  table  wherein  range  of  the  objective  function  F,-  is 
determined  by  the  minimum  and  the  maximum  values  for  the 


e-constraint  technique.  So,  the  range  of  each  objective  function  is 
determined  as  Eq.  (2): 

r,  =  F“ax  —  F“in  (2) 


The  next  step  after  determining  the  range  (r,)  pertaining  to  all 
objective  functions  from  payoff  table  is  to  divide  the  range  of  the 
remaining  objective  functions  F2,...,  Fp  to  n2,...,  np  into  same  intervals 
utilizing  (n2  —  1 ),. . .,  (np-l)  intermediate  grid  points  with  equal 
distances,  respectively.  The  total  grid  points  with  consideration  of 
the  minimum  and  the  maximum  values  of  the  range  would  be 
(n2+l ),...,  (np+l)  corresponding  to  F2,...,  Fp,  respectively.  Therefore, 
the  Pareto  optimal  solution  is  derived  through  solving  (n2+l)x 
(n3+l)x  ...  x  (np+l)  optimization  sub-problems.  The  importance  of 
objective  functions  is  excluded  by  the  augmented  e-constraint  tech¬ 
nique  once  the  Pareto  solutions  are  generated  which  is  incompatible 
with  the  DM  purposes.  However,  augmented-weighted  e-constraint 
technique  is  used  in  this  paper  so  that  the  relative  significance  of  the 
objective  functions  is  clearly  modeled  in  the  Pareto  solutions  genera¬ 
tion.  In  this  regard,  the  objective  function  in  (1)  is  changed  as  below 
according  to  the  augmented-weighted  e-constraint  technique: 


Min  /  Max 


WiF^xj+dir,  r, 


s.t.  Fi(x)-diri  Sj  =  e,  V  s,  e  R+ ,  i  =  2,3,...,p.  (3) 

where,  the  weighting  factors  of  DM  for  the  ith  objective  function 
are  denoted  by  w,-.  Note  that,  the  technique  presented  here  as  augm¬ 
ented-weighted  e-constraint  is  quite  different  from  the  procedure  of 


Fig.  1.  Flowchart  of  the  e-constraint  optimization  method  for  the  MMP  problem. 
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weighting  method  that  results  in  a  single  solution  for  a  set  of  weights 
without  ensuring  the  efficiency  of  the  solution.  The  proposed  math¬ 
ematical  formulation  in  (4)  must  be  solved  to  generate  Pareto  optimal 
solution  in  each  sub-problem. 


Min  /Max 


F  i  (x)  + 


dir,  i-i 

Wi 


£  WjS? 
i  =  2  ri 


s.t.  Fi(x)-diri  sf  =  ef  V  s,-  e  R+,  i  =  2,3 ,...,p,k  =  0, 1. 


(4) 


„k  cMin (dirt  + 1 )  cMax(diri  - 1 )  ,  (din  n  k) 

'  '  2  2  +  (jj 

V  i  =  2,3,  ...,p,  l<  =  0, 1....,  q,-  (5) 


2.2.  Fuzzy  decision  maker 


A  Fuzzy  decision  making  technique  is  used  in  this  paper  as  Refs. 
[34,35]  after  deriving  the  Pareto-optimal  solutions  through  solving 
the  optimization  sub-problems  to  select  the  most  preferred  solu¬ 
tion  in  accord  with  particular  purpose  of  the  application.  In  this 
technique,  for  each  objective  function  of  the  MMP  problem,  a 
linear  membership  function  is  defined.  Eq.  (6)  is  defined  as  below 
for  objective  functions  that  are  set  to  be  minimized: 


r  i 


Mi  =  • 


0 


rr  rmin 

J  i  —  Ji 

rr  rmax 

Ji  —Ji 


(6) 


And  for  the  objective  function  intended  to  be  maximized,  the 
linear  membership  function  is  defined  in  (7),  as  well. 


(  0 


Mi  =  ' 


fi  -fTm 

/r-/r 

i 


/■  </rin 

rr  <fi  </r 

rr  rmax 

Ji  — Ji 


(7) 


where,  the  range  of  objective  function/,  derived  from  payoff  table 
is  denoted  by/™n  and  / jTlax.  The  value  of  the  ith  objective  function 
in  the  rth  Pareto  optimal  solution  as  well  as  corresponding 
membership  function  are  indicated  by  /[  and  p\,  respectively. 
The  concept  pursued  by  the  value  of  obtained  membership 
function  is  the  nicety  of  the  solution  derived  for  the  ith  objective 
function  in  the  rth  Pareto  optimal  solution.  Besides,  for  the  rth 
Pareto  optimal  solution,  the  total  membership  function  (/)  is 
defined  on  the  basis  of  its  individual  membership  functions  as  (8): 


M  = 


(8) 


where,  ith  objective  function's  weighting  in  the  MMP  problem  is 
denoted  by  w,  and  the  p  indicates  the  number  of  objective 
functions.  By  considering  the  significance  of  economic  and  secur¬ 
ity  aspects,  the  system  operator  can  choose  the  weighting  value  w,. 
The  best  Pareto  optimal  solution  of  the  MMP  problem  is  the  one 
with  the  highest  value  of  membership  function  //  with  respect  to 
the  assigned  weight  factors. 

The  presented  MMP  solution  method  is  formed  by  coming 
together  the  augmented-weighted  e-constraint  technique  and  lex¬ 
icographic  optimization.  The  procedure  of  the  proposed  method  can 
be  stated  as  follows: 


Step  1 :  By  employing  the  lexicographic  optimization  approach, 
the  payoff  table  pertaining  to  a  MMP  problem  is  computed. 
Step  2:  The  range  of  the  ith  objective  function  (i= 2,  3,...,  p )  is 
determined  using  payoff  table  as  Eq.  (2). 


Step  3:  According  to  the  formulation  proposed  in  (5),  the  range 
of  at  least  p-1  objective  functions  is  divided  into  q,  (i= 2,  3  ,...,p) 
equal  intervals. 

Step  4:  The  feasible  optimization  sub-problems  in  (4)  are 
solved  applying  the  presented  MMP  solution  method  to  pro¬ 
duce  the  Pareto  efficient  solution  while  the  infeasible  ones  are 
discarded. 

Step  5:  The  efficient  solutions  derived  from  step  4  are  evaluated 
using  the  Fuzzy  decision  making  process  (6)-(8)  to  choose  the 
most  desired  Pareto  optimal  solution. 

The  solution  method  on  the  basis  of  e-constraint  optimization 
technique  for  the  proposed  MMP  problem  is  depicted  in  Fig.  1 
as  a  flowchart. 


3.  MIP  formulation  for  SCUC 


The  two  objective  functions  of  the  presented  framework  for  the 
SCUC  problem  can  be  stated  as  (9): 


MultiObjective  function 


p!  cos  t  min  imization 
F2  emmision  min  imization 


(9) 


where,  the  detailed  description  of  the  objective  functions  F ),  and 
F2  pertaining  to  the  SCUC  problem  is  given  in  the  following 
section. 


3.2.  Cost  minimization 

Minimizing  cost  is  the  main  objective  function  in  the  proposed 
SCUC  problem  as  (10): 

Ft  :MinFC0ST=  £  (  T,Ajy(j,  t)+  XAZ(i,  t)+B(i,  t)  +  F(i,  t)  +  C(i,  t) 

t  eT\J  ej  i  el 

(10) 

where,  the  total  cost  is  denoted  by  FC0ST  including  thermal  and 
hydro  units  operation  cost  that  it  is  intended  to  be  minimized  by 
ISO.  The  start-up  cost  is  indicated  by  the  first  term  [8,36]  while 
thermal  unit’s  cost  comprising  fuel  cost,  valve  loading  effects  cost, 
shut-down  cost  and  start-up  cost  will  be  explained  in  detail  later. 
Because  of  wear  and  tear  of  the  windings  and  mechanical  equip¬ 
ment  as  well  as  loss  of  water  during  start-up  and  maintenance  in 
addition  to  malfunctions  in  the  control  system,  start-up  cost  of 
hydro  units  is  considered  by  authors  in  Ref.  [36]. 


3.2.  Emission  minimization 


F2  :  Min  FEMISS,0N 


£  £ 

i  e  EMGt  e  T 


emin(0 Id(U  t)+SUE(i)y(i,  t)+SDE(i)z(i,  t ) 

NE 

+  £  \Pn(F t)E(p“_,(i))+be„(i)d„(i,t)] 

n  =  1 


(if) 


where,  the  first  term  indicates  the  emission  caused  by  thermal 
units  which  are  off  and  provide  non-spinning  reserve  [37]  and 
EMC  =  [S02,  NO*)  as  the  most  significant  emissions  caused  by 
electricity  generation  section  affecting  the  environment  are  S02 
and  NO*  [28,38].  Authors  of  Ref.  [28]  have  modeled  the  emissions 
as  a  quadratic  term  while  the  linearized  form  using  piecewise 
linear  approximation  is  illustrated  in  Fig.  2.  It  should  be  noted  that, 
in  order  to  more  precise  modeling,  emission  caused  by  the  start¬ 
up  and  shut-down  of  thermal  units  is  included  in  the  emission 
function. 
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Emission  Generation  (lbs/h) 


Fig.  2.  Piecewise  linear  emission  generation  curve  with  M  prohibited  zones. 


Generation  ( MW ) 


Fig.  3.  Piecewise  linear  non-concave  unit  performance  curve  for  hydro  unit  j  at 
hour  t. 


3.3.  Hydro  units'  model 

The  relationship  between  water  discharge,  generated  power 
and  multi-head  of  reservoir  is  depicted  in  Fig.  3  where  the 
connection  of  hydro  units  is  both  parallel  and  in  series.  Multi- 
head  reservoir  is  used  due  to  the  small  storage  capacity  of 
reservoirs  and  also  because  the  generated  power  depends  on  the 
hydro  unit  head.  The  number  of  heads  is  shown  with  L  in  the 
following  formulations  while  the  higher  number  of  L  improves  the 
model  accuracy. 

3.3.  J.  Linear  formulations  for  volume  and  multi-head 

This  section  comprises  the  linear  formulations  of  hydro  power 
units  with  L  heads  (Fig.  3): 

>  ■■■Vj  eJ,Vt  e  T  (12) 

V(j,  t)  <  VL(j)PL- 1  O’,  t)  +  £  V„_,O1|)0„-2O'.t)-A,-lO'.t)]  Vjej,  VteT 

n  —  2 

(13) 

v(j,t)>  Vt_i<M-l<J,t)+  X  Vjej, VteT 

n  —  3 

(14) 

(15) 


where, 

A>(7.t)  =  l. 

when  n  +  lth  is  used,  binary  variable  pn(j,t)  would  be  equal  to  1. 
Constraint  (12)  eliminates  the  0-1  conflicts  between  these  binary 
variables.  Right  head  corresponding  to  volume  value  is  indicated 
by  constraints  (13)-(15).  Constraints  (13)-(15)  ensure  that  the 
volume  value  of  each  hydro  unit  at  each  time  is  greater  than  the 
minimum  content  of  that  hydro  unit. 

3.3.2.  Linear  power-discharge  curve 

The  linear  relationship  between  discharged  water,  generated 
power  as  well  as  variable  head  is  declared  as  below: 

P(j-t)-pk(j)Kj,t)-  £  £  A.0,1) 

neN  n = 1 

+  LZ  0]  <  o  Vj  ej,  Vt  e  T,  1  <  k  <  L  (16) 

n  =  k 

P(i,t)-p.<j)I<J,t)-  £  q„(j,t)  bkn(j)+p(j)[(l<-1)-  £  Pn(j,  t) 

neN  n  =  1 

+  1£  find,  t)]>OVjeJ,VteTA<k<L  (17) 

n  =  k 

The  power  generated  by  hydro  unit  j  at  time  t  and  the 
minimum  power  generation  of  head  k  are  denoted  by  PQ'.t)  and 
pfc0).  respectively.  pn(j,t )  indicates  the  proper  head  and  p(j)  is  the 
capacity  of  hydro  unit  j.  qn(j,  t )  denotes  the  discharge  of  block  n 
while  the  slope  of  the  block  n  pertaining  to  the  performance  curve 
k  of  hydro  unit  j  is  represented  by  b„(j). 

3.3.3.  Water  discharge  limits 

Q(j,  t)  =  am,  ty +  £  q„(j,  t)  Vj  ej,  V  t  e  T  (18) 

neN 

9(j,  t)<(f(j,  t)  +  S(j,  t)  <  0(j,  t)  VjeJ.VteT  (19) 

where,  water  discharge  of  hydro  unit  j  at  time  t  is  indicated  by 
Q(j,t )  and  equals  to  the  minimum  water  discharge  provided  that 
the  hydro  unit  is  on,  plus  the  water  discharge  of  blocks.  However, 
in  order  to  meet  flood  and  irrigation  necessities,  water  discharge 
must  be  limited. 

s(j,t)=  £  bn(j)yfn(j,t)  Vjej,  VteT  (20) 

n  =  1 

v(j,t)  =  v0(j)+  2  yr„(j,t)  VjeJ.WteT  (21) 

n  =  1 

[Vi  (J)  -  vo(j)]Pi  Q,  t)  ^  ¥\  (j,  t)  <  Vi  0j  -  VoOj  Vjej.VteT  (22) 

[Vn(j)-v„-i(j)Wn(j,t)<¥n(j,t)  n  =  2,3,....,L,VjeJ,VteT  (23) 

Wnij,  t)  <  [v„0j  -  vn_  t  (j)]Pn  _  i  O',  t)  n  =  2, 3, ....,  L,  Vjej,  VteT  (24) 

The  reservoir  spillage  is  modeled  as  a  function  of  reservoir 
volume  value.  The  slope  of  the  block  n  pertaining  to  the  reservoir 
volume  value  is  denoted  by  bn(j)  and  block  n  of  the  reservoir 
volume  value  is  denoted  by  y/„(j,  t).  Constraint  (21 )  indicates  the 
relationship  between  reservoir  content  and  the  volume  blocks.  The 
first  block  of  the  reservoir  volume  is  restricted  by  constraint  (22) 
while  other  reservoir  volume  blocks  at  each  time  are  restricted 
through  constraints  (26)-(28). 

q  1 0,  t)  <  QiOVO'.t)  vjej,  VteT  (25) 


v(j,  t)  >  v0(j)  Vjej,  VteT 
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giO'.t)>Qi(/)friO'.t)  Vj  e-/j  VteT  (26) 

qnO\t)<Qn(j)h„-,(j,  t)  Vjej, VteT, VneN  (27) 

Q„<j,  £)  >  Q„(j)hn(j,  t)  Vj  ej,  VteT,  VneN  (28) 

Q.„0)  indicates  the  maximum  water  discharge  of  the  nth  block 
pertaining  to  hydro  unit  j  and  the  binary  variable  hn(j,  t )  is  equal  to 
1  if  the  water  discharge  of  hydro  unit  j  at  time  t  has  overstepped 
the  nth  block.  Constraints  (25)  and  (26)  restrict  the  water 
discharge  of  the  first  block  while  other  blocks  of  water  discharge 
are  restricted  through  constraints  (27)  and  (28). 

3.3.4.  Initial  and  final  volume 


v(j,0)  =  v°(j) 

VjeJ 

(29) 

v(j,  0)  =  vs(j) 

Vjej 

(30) 

The  initial  and  the  final  volume  of  the  reservoir  at  the  begin¬ 
ning  and  the  end  of  the  planning  horizon  are  limited  through 
constraints  (29)  and  (30). 

3.3.5.  Water  balance 

In  the  cascaded  configuration  of  hydro  units,  the  water  dis¬ 
charged  from  the  upstream  hydro  unit  is  the  inflow  of  the  next 
downstream  hydro  unit.  In  Eq.  (31)  which  indicates  the  continuity 
equation,  the  forecasted  natural  water  inflow  to  the  reservoir  of 
the  hydro  unit  j  at  time  t  is  indicated  by  F(j,t)  and  ;/  denotes  the 
conversion  factor  which  is  utilized  to  convert  water  discharge 
reservoir  into  stored  water.  Qj  denotes  the  upstream  reservoirs  of 
hydro  unit  j  and  the  time  delay  corresponding  to  the  water 
transport  from  hydro  units  i  to  j  is  expressed  by  r,j. 

v(j,  t)  =  v(j,  t  —  1 ) + F(j,  t)  —  t/.(Q(j,  t)+s(j,  t)) 

+  4  I  [Q(/,t-ry)  +  S(/,t-Ty)]  VjeJ.VteT  (31) 

i  e  £2j 


3.3.6.  Operating  services 

Since,  the  hydro  units  are  known  as  fast  response  units,  they 
would  be  suitable  to  provide  spinning  reserve  as  well  as  non¬ 
spinning  reserve.  This  paper  uses  formulations  proposed  in  Refs. 
[37,39]  as  below: 

p  (j,  t)  +  R(j,  t)  +  Nu(j,  t)<p  (j)  l(j.  t)  VjeJ.VteT  (32 ) 

Nd(j,  t)  <  QSC(j)(  1  -  I(j,  t))  VjeJ.VteT 

3.3.7.  Logical  status  of  commitment 

The  logical  relationship  between  three  binary  variables  corre¬ 
sponding  to  start-up  status,  shut-down  status  and  the  status  of 
unit  at  each  hour  of  the  planning  horizon  is  represented  in 
constraint  (33)  while  constraint  (34)  eliminates  the  conflicts 
between  simultaneous  start-up  and  shut-down  by  the  unit. 

y(j,  t)—Z(j,t)  =  I(j,t)-I(j,t-1)  VjeJ.teT  (33) 

y(j,  t)+Z(j,  t)  <  1  VjeJ.teT  (34) 

3.4.  Thermal  units'  model 

The  linear  model  of  thermal  units  is  proposed  through  below 
formulations: 


3.4.1.  Fuel  cost  function  considering  POZ 

Generally,  fuel  cost  function  of  thermal  units  is  presented  by  a 
quadratic  function  while  the  fuel  cost  function  pertaining  to 
thermal  units  is  discrete  as  they  are  not  able  to  operate  in  specific 
zone  due  to  physical  operating  restriction.  The  equivalent  piece- 
wise  linear  model  of  the  fuel  cost  function  is  depicted  in  Fig.  4 
with  MPOZs. 

M  +  l 

F(i,t)=  2  [/in(i,  t)F(p"_ ,(!))  +  b„{i)8n(i,t)]  Vi  el,  Vt  eT  (35) 

n  =  1 

M  +  l 

P(i,t)=  I  [P!i-1(i)/5n(i.t)  +  5n(i.t)]  Viel,  VteT  (36) 

n  =  1 

F(i,t )  is  piecewise  linearized  fuel  cost  function  and  binary 
variable  p„(i,  t)  is  equal  to  1  if  power  block  n  for  unit  i  of  piecewise 
linearized  fuel  cost  curve  is  selected.  Cost  of  generation  p“  _-,(!)  is 
equal  to  F(p“  1  (i)).  Slope  and  generation  of  power  block  n  are  used 
at  second  term.  The  power  produced  by  thermal  unit  i  at  hour  t  is 
the  sum  of  the  “power  length”.  Other  constraints  for  linearization 
of  fuel  cost  function  can  be  stated  as  follows  [40]: 

The  piecewise  linearized  form  of  the  fuel  cost  function  is  denoted 
by  F(i,t )  while  JSJi,  t )  is  a  binary  variable  which  is  equal  to  1  provided 
that  the  power  block  n  for  unit  i  of  piecewise  linearized  fuel  cost 
curve  is  selected.  Moreover,  generation  cost  denoted  by  p“  _  1  (i)  is 
equal  to  F(p“  1  (z)).  The  second  term  comprises  slope  and  genera¬ 
tion  pertaining  to  power  block  n.  the  power  generated  by  thermal 
unit  i  at  time  t  is  the  sum  of  the  “power  length”.  Other  constraints 
considered  for  linearizing  the  fuel  cost  function  are  stated  as 
below  [40]: 

Sn(i,t)>0  n=  1,2,..., M+l, Viel,  VteT  (37) 

8n(i,t)<\pd„(.i)-p“-i(j)Wn(i,t)  n=  1,2,..., M  +  l,  Viel,  VteT 

(38) 

M+l 

I  Pn(i,t)  =  l(i,t)  Viel,  VteT  (39) 

n  =  1 

P„(i,t)e{ 0,1)  n  =  l,2 . M+l, Viel,  VteT  (40) 

where,  p%(i)  =  pmin(i)  and  p£,+1(!)  =pmax(i). 

The  power  generated  at  each  block  should  be  positive  and 
limited  to  its  upper  bound.  The  operation  of  the  thermal  units  at 
only  one  of  the  operating  zones  is  implied  in  constraint  (39)  e.g.  if 
block  2  of  the  fuel  cost  curve  is  chosen  for  the  thermal  unit 
i  (/?2(i,  t)  is  equal  to  1),  afterward  S2(i,t)  can  be  less  or  equal  to  the 
maximum  “power  length”  of  the  second  block  (p^(i) — p“(i))  and 
other  "other  lengths”  are  set  to  zero.  Consequently,  thermal  unit 
generation  at  this  time  would  be  p'j(i)  plus  82(i,t). 


Fuel  cost  (Mbtu/h) 
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3.4.2.  Valve  loading  effects'  cost 

In  thermal  generating  units  with  high  installed  capacity,  there 
are  multi-valve  steam  turbines.  Once,  the  generated  power  is 
needed  to  be  increased  in  such  thermal  units  by  ISO,  the  multi¬ 
valve  turbines  are  opened  in  order;  As  a  result,  a  ripple-like  impact 
is  produced  in  the  heat  rate  curve  of  the  thermal  units  [41],  The 
valve  loading  effects  are  taken  into  consideration  in  Refs.  [41,42]  as 
an  absolute  sinus  function  of  the  generated  power  while  applying 
this  sinus  function  to  the  model  brings  severe  challenges  in  solving 
the  realistic  and  large  power  systems  due  to  high  degree  of  non¬ 
convexity  and  also  the  state  of  being  non-smooth.  Thus,  linear 
formulation  is  proposed  in  this  paper  to  deal  with  the  valve 
loading  effects  (Fig.  5). 


C(i,  t)  = 


2e/j 


s/2  t  [^n  +  lO'.tl-V^n+^ff)] 


+  (2-\/2)  2  Wun+zlUt) -yzAn+3(i,t)\ 
n  =  0 


Vie/,  VteT 


(41) 


P(i,t)  =  pmin(i)I(i,  t)+ 


ki 

Z 


n  =  0 


V/4n  +  l(b  0  +  V'4n  +  2(,»  0 

+  W4n+3(U  t)  +  l//4„  +  4(i,  t) 


Viel,  VteT 
(42) 


P(i,  t)  <  pmax(i){I(i,  t)  -Z(i,t+ 1 )}  +  SD(i)Z(i,  t+ 1 )  Vie  /.VteT 

(47) 

p(i,t-1)-p(i,t)<SD(i)Z(i,t)  +  RDL(p(i,t))  Viel, VteT  (48) 


p(i,  t+ 1  )-p(i,  t)  <  SU(i)y(i,t+L)+RUL(p(i,  t))  Viel, VteT 

(49) 

The  power  output  limits  of  thermal  units  are  given  in  con¬ 
straint  (46)  implying  that  the  generated  power  at  time  t  should  be 
less  than  the  upper  bound  and  greater  than  the  minimum  output 
power  provided  that  the  unit  is  committed.  The  upper  bound  of 
thermal  units  at  any  time  of  study  is  indicated  in  constraint  (47) 
considers  the  status  of  thermal  units.  If  the  status  of  the  unit  is  off 
at  next  hour,  the  generated  power  should  be  less  than  the  shut¬ 
down  limit.  The  shut-down  ramp  rate  as  well  as  Ramp-Down  Limit 
(RDL)  is  represented  in  constraint  (48).  Constraint  (49)  includes 
the  start-up  ramp  rate  and  Ramp-Up  Limit  (RUL).  Furthermore,  the 
dynamic  ramp  rate  is  considered  in  this  paper  instead  of  fixed  rate 
to  derive  more  precise  model  unlike  majority  of  published  papers. 


4*ri(i,t)<VU(L  t)<4^./(i,  t)  Vi  el, VteT  (43) 

t)  <  i/r„(i,  t)  <  ^Xn  - 1 0',  0  Viel, VteT, n  =  2.3,...,x,  (44) 

Xn(i,t)e{  0,1}  Viel, VteT, n  =  \,2,...,Xi  (45) 

where,  k,-  =  floor[fi(pmax(i)~pmin(i)A)]  and  x,  =  floor[4/i(pmax(i)  - 

Pinin  *  9  zc)  ]. 

The  valve  loading  effects’  cost  is  denoted  by  C(i,t )  and  the 
coefficients  of  valve  point  effects  pertaining  to  the  thermal  unit  i  at 
time  t  are  represented  by  f  and  e,  while  i//  n(i,  t )  indicates  the 
power  generated  by  block  n.  The  thermal  power  produced  by 
thermal  unit  i  at  time  t  would  be  the  sum  of  minimum  output 
power  once  the  unit  is  committed  plus  power  generated  in  each 
block.  The  power  generated  in  the  first  block  is  limited  by 
constraint  (43).  This  generated  power  must  be  smaller  than  or 
equal  to  nl4f  and  greater  than  or  equal  to  zero  i.e.  “power  length” 
of  each  block.  /( i, t)  ensures  that  if  unit  i  is  decommitted  at  time  t, 
then  the  power  produced  at  the  first  block  is  zero.  Constraints  (43) 
and  (44)  are  used /„(!,  t )  to  limit  the  power  generated  in  each  block 
introduce  binary  variable.  Each  time  the  power  generated  by  the 
thermal  unit  i  at  time  t  has  overstepped  block  n,  this  binary 
variable  is  equal  to  1. 


3.4.4.  Dynamic  RDL  and  RUL 

The  ramp  rate  limits  restrict  the  power  generated  between  two 
successive  operating  periods.  The  generators  respond  to  the  hourly 
variations  in  the  system  demand  through  increasing  or  decreasing 
the  produced  power.  The  presented  dynamic  ramp  rate  is  a 
function  of  the  output  power  taking  into  consideration  POZs 
(Fig.  6). 


M+l 


RDL(p(i,t))=  2  RDL„(i)Pn(i,t) 

n  =  1 

Vie  I,  VteT 

(50) 

M  +  l 

RUL(p(i,t))=  2  RULn(i)/]n(i,t) 

Vi  el,  VteT 

(51) 

n  =  1 


The  dynamic  RDL  and  RUL  are  represented  in  constraints  (50) 
and  (51),  respectively  while  the  selected  operating  zone  is  denoted 
by  //„(i,  t) . 

3.4.5.  Time  varying  start-up  cost  function 

The  cold  start-up  cost  can  be  represented  by  the  exponential 
function  of  number  of  hours  that  a  thermal  unit  is  off.  Ref.  [7]  uses 
linear  form  of  varying  start-up  cost  function. 


3.4.3.  Generation  thermal  unit  capacity  limits 

Pmin(i)Ri,  t)  <  p(i,  t)  <  p(i,  t)  Vi  el, VteT  (46) 


Valve  loading  cost  ($) 


3.4.6.  Reserve  services 

Refs.  [43,44]  use  operating  services  to  support  sudden  events 
like  outages  of  transmission  lines  and  generators.  Operating 
reserves  include  spinning  and  non-spinning  reserves.  Spinning 
reserve  is  defined  as  the  unloaded  units  of  synchronized  genera¬ 
tion  being  available  for  use  in  ten  minutes  while  non-spinning 
reserve  is  defined  as  the  unsynchronized  generating  units  that  are 
available  for  use  in  ten  minutes.  The  formulations  used  in  this 
paper  for  spinning  and  non-spinning  reserve  are  the  same  as  the 
ones  proposed  in  Ref.  [44]. 

3.4.7.  Minimum  up-time  (MUT)  and  minimum  down-time  (MDT) 

If  a  unit  is  turned  on,  it  should  not  be  turned  off  for  at  least 
several  hours,  as  it  is  restricted  by  the  MUT  limitation.  On  the 
other  hand,  if  a  unit  is  turned  off,  it  should  not  be  turned  on  for  at 
least  several  hours.  In  this  regard,  formulations  presented  in  Refs. 
[7,39]  are  used  in  this  paper. 
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3.4.8.  Logical  status  of  commitment 

The  formulation  proposed  in  section  2.3.7  [7,39]  is  utilized  here 
to  eliminate  the  conflicts  between  start-up  and  shut-down  status 
at  the  same  time. 

3.5.  Network  security 

The  transmission  constraints  are  modeled  as  linear  constraints 
on  the  basis  of  DC  power  flow  model  as  [45], 


found  in  Ref.  [8]  in  which  thermal  units  5, 10, 11,  28,  36,  43,  44  and 
45  have  valve  loading  effects  cost  and  thermal  units  7, 10,  30,  34, 
35  and  47  have  POZs  limitations.  Thermal  units  8,  9, 10, 11, 15, 16, 
17,  18,  21,  22,  23,  24,  28,  29,  30,  31,  32,  33,  34,  35  and  36  can 
generate  emission.  The  data  pertaining  to  the  demand  and  thermal 
units  is  available  in  Ref.  [46],  The  computer  system  employed  in 
this  paper  to  solve  the  M1P  optimization  of  SCUC  problem  is  a 
Pentium  IV  PC  with  2.6  GHz  clock  speed  and  2  GB  RAM  using 
solver  CPLEX  in  the  GAMS  [47]  environment. 


4.  Case  study 


4.1.  Numerical  results 


The  model  is  implemented  on  the  modified  IEEE  118-bus 
system  as  shown  in  Fig.  7  to  verify  the  effectiveness  of  the 
proposed  framework.  This  system  comprises  54  thermal  units 
(33  coal-fired  units,  11  gas-fired  units  and  10  oil-fired  units)  as 
well  as  8  hydro  plants.  The  detailed  data  of  this  system  can  be 


Rampuplimit(MW) 


RULm.i 


(O 


In  this  paper,  lexicographic  optimization  and  hybrid  augmented 
e-constrain  are  used  to  generate  the  Pareto  solutions  of  SCUC.  The 
objective  function  Fj  i.e.  cost  minimization  is  taken  into  considera¬ 
tion  as  the  main  objective  function  in  the  E-constraint  technique. 
49  grid  points  (q2=49)  is  selected  for  F2  (objective  function 
permssion-j  t0  derive  Pareto  optimal  solutions  from  the  MMP  pro¬ 
blem.  The  Pareto  solutions  are  generated  by  solving  the  problem 
for  50  times  (q2+l=50)  while  all  solutions  are  feasible  [15,27,35], 
Two  cases  are  considered  to  study  the  presented  SCUC  where  the 
first  case  only  uses  augmented  E-constraint  technique  without 
utilizing  lexicographic  optimization.  The  resulted  payoff  table  (Cq ) 
is  represented  as  below: 


4>i  = 


172,224.770 

348,006.792 


54,417.601 

1908.816 


Pmin(i)  Pd‘(i)  p„'(i)  Pd2(i)  P„M(i)  p„.,(i) 

Generation(MW) 

Fig.  6.  Ramping  up  limit  as  a  stepwise  linear  function  of  generation  with  M  POZs. 


The  first  column  corresponds  to  the  first  objective  function 
(Fi=Fcosr)  and  the  second  column  gives  the  results  obtained  for 
the  objective  function  F2  i.e.  emission  minimization  while  all  the 


System  Description: 

118  buses 
186  branches 
91  load  sides 
54  thermal  units 


One-line  Diagram  of  IEEE  118-bus  Test  System 
I IT  Power  Group,  2003 


Fig.  7.  Single-line  diagram  of  118-bus  system. 
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obtained  results  are  in  terms  of  $.  Moreover,  the  main  diagonal  of 
the  payoff  table  <£]  represents  the  results  obtained  from  individual 
optimization  of  the  objective  functions  (the  Utopia  points).  As  it 
can  be  observed,  the  single  objective  SCUC  with  F i  as  the  objective 
function  leads  to  the  minimum  cost  but  high  amount  of  emission 
is  generated  in  this  case  (54417.601  lbs).  On  the  other  hand, 
the  results  obtained  for  individual  optimization  of  F2  shows  much 
less  emission  generation  (1908.816  lbs)  but  with  higher  cost 
(348006.792  $). 

r  172224.770  36642.31' 

02  =  232511.970  1908.816 


As  it  can  be  seen  in  the  values  of  the  first  row  of  payoff  tables 
0!  and  0 2 .  the  Nadir  point  pertaining  to  the  second  objective 
function  F2  is  enhanced  by  17775.291  (lbs)  i.e.  (54417.601- 
36642.31  =  1775.291  (lbs)).  Therefore,  narrower  ranges  are  derived 
for  all  objective  functions  in  02  compared  to  0 ^ ,  as  the  second  row 
of  the  02  includes  better  results  than  those  in  0\.  This  is  because 
some  inefficient  Pareto  solutions  generated  using  e-constraint 
technique,  are  discarded  by  lexicographic  optimization  and  hybrid 
augmented  e-constraint  method. 

It  is  implied  from  payoff  table  02  that,  cost  can  be  decreased  by 
ISO  to  172224.770  ($),  if  discards  emission  generation  but  results  in 
emission  equals  to  36642.31  (lbs).  On  the  other  hand,  if  the  second 
objective  function  is  only  taken  into  consideration  i.e.  pem,sswni  the 
amount  of  emission  generation  would  decrease  to  1908.816  (lbs)  by 
ISO  but  leading  to  increase  in  cost  to  232511.970  ($).  If,  ISO  have  to 
consider  the  minimization  of  the  emission  generation  in  its  multi¬ 
objective  SCUC  problem,  the  cost  would  increase  by  60287.2  ($)  i.e. 
(232511.97-172224.77=60287.2  ($))  in  spite  of  substantial  reduc¬ 
tion  in  emission  generation  by  34733.494  (lbs)  i.e.  (36642.31- 
1908.816=  34733.494  (lbs)).  This  state  shows  the  conflicting  nature 
of  the  two  objective  functions  proposed  in  this  paper. 

A  Fuzzy  decision  maker  is  utilized  in  this  paper  as  Ref.  [34]  to 
select  the  most  compromise  solution  among  the  Pareto  solutions 
derived  for  the  multi-objective  SCUC  problem.  In  the  proposed 
Fuzzy  decision  maker,  weighting  factors  corresponding  to  two 
objective  functions  are  taken  as  the  same.  The  total  membership  of 
all  Pareto  solutions  is  depicted  in  Fig.  8  while  Pareto  solution  40  is 
selected  as  the  best  one  due  to  higher  value  of  total  membership 
obtained  for  this  Pareto  solution  as  0.663.  Table  1  represents  the 
results  obtained  for  the  Pareto  solution  number  40  considering 
equal  weighting  factors.  The  total  membership  value  in  Table  1 
shows  the  degree  of  optimality  and  equals  to  0.663.  The  most 
preferred  solution  can  be  easily  found  by  the  DM  by  changing  the 
weighting  factors  while  ISO  explicitly  tends  to  minimize  the  cost 
rather  than  the  emission.  Therefore,  by  selecting  the  weighting 
factors  as  2  and  1  for  cost  and  emission,  respectively,  the  obtained 
solutions  in  this  state  are  represented  in  Table  2  while  it  is  obvious 
that  a  Pareto  solution  with  higher  value  of  cost  membership  and 


Fig.  8.  Total  membership  value  of  Pareto  solutions  for  equal  weighting  factor  case. 


Table  1 

Optimum  solution  of  SCUC  with  equal  weighting  factor. 


Objective 

function 

Weighting 

factor 

Objective  function 
value 

Membership 

value 

Cost 

1 

200520.56 

0.772 

Emission 

1 

8997.28 

0.796 

Total  membership  of  all  objective 

0.663 

functions 

Table  2 

Optimum  solution  of  SCUC  with  different  weighting  factor. 

Objective 

Weighting 

Objective  function 

Membership 

function 

factor 

value 

value 

Cost 

2 

175381.22 

0.948 

Emission 

1 

31680.38 

0.143 

Total  membership  of  all  objective 

0.679 

functions 

no] 

_ 

20  30 

Pareto  solution  number 

“Cost  membership  1  Emission  membership  — 


Total  membership 


Fig.  9.  Variation  of  objective  functions  and  total  membership  value  versus  Pareto- 
optimal  solutions  for  different  weighting  factor  case. 


Table  3 

Optimization  statistics  for  multi-objective  SCUC  problem. 


No.  of  single 

No.  of  single 

No.  of  discrete 

No.  of 

Solution 

constraints 

variables 

variables 

iteration 

time  (s) 

1521  666 

1309  986 

566  676 

512  432 

986 

1200 


1  2  3  4  5  6  7  8  9  10  11  12  13  14  15  16  17  18  19  20  21  22  23  24 

Time  (h) 

■  Phlt  «Ph2t  nPh3t  nPh4t  aPsIt  aPs2t  sPs3t 


Fig.  10.  The  details  of  the  best  compromise  solution  with  equal  weighting  factor. 

lower  value  of  emission  membership  will  be  searched.  According 
to  Table  2  and  Fig.  9,  a  remarkable  improvement  in  cost  member¬ 
ship  is  observed  as  0.772  in  the  case  of  equal  weighting  factors  to 
0.948  in  the  case  of  different  weights.  Table  2  represents  the  cost 
value  as  175381.22  ($)  which  is  very  close  to  its  ideal  value 
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Table  4 

Results  obtained  from  different  methods. 


Method 

Ref.  [20] 

Ref.  [29] 

Ref.  [30] 

Ref.  [31] 

Proposed  method 

ELS 

EES 

CEES 

ELS 

EES 

CEES 

ELS 

EES 

CEES 

ELS 

EES 

CEES 

ELS 

EES 

CEES 

Cost  (S) 

NR 

NR 

47,906 

43,500 

51,449 

44,914 

41,909 

45,392 

43,507 

43,278 

47,871 

44,344 

40,766.83 

41,145.14 

40,854.28 

Emission  (lbs) 

NR 

NR 

26,234 

21,092 

18,257 

19,615 

30,724 

17,659 

18,183 

17,984 

17,019 

17,408 

18,278.57 

15,666.61 

16,180.20 

Solution  time  (s) 

NR 

NR 

4,582 

72.96 

72.74 

74.97 

NR 

NR 

NR 

NR 

NR 

NR 

1.78 

1.12 

15.68 

NR=not  reported. 

i.e.  172224.77  but  higher  value  of  emission  is  obtained  i.e. 
31680.38  (lbs)  while  its  membership  value  is  low  (0.143).  The 
variation  of  objective  functions  and  total  membership  value  vs. 
Pareto-optimal  solutions  corresponding  to  different  weighting 
factor  case  is  shown  in  Fig.  9. 

4.2.  Discussion 

An  efficient  multi-objective  method  is  proposed  in  this  paper 
that  can  be  used  by  ISO  for  day-ahead  scheduling  of  generation 
units  even  for  large  systems  with  lots  of  variables  and  constraints. 
The  detailed  optimization  statistics  for  the  multi-objective  SCUC 
problem  solved  by  lexicographic  optimization  and  hybrid  aug¬ 
mented  e-constraint  method,  is  represented  in  Table  3.  So  far, 
there  is  no  published  paper  implementing  the  SCUC  problem  for 
such  a  large  system  (IEEE  118-bus)  comprising  54  thermal  units  as 
well  as  8  hydro  plants  while  two  conflicting  objective  functions  as 
cost  and  emission  are  optimized,  simultaneously.  In  addition, 
many  practical  constraints,  such  as  network  security,  valve  loading 
effects,  dynamic  ramp  rate  limit  and  multi  POZs  of  thermal  units 
and  multi  performance  curves  for  hydro  units  are  considered  in 
this  paper.  As  a  result,  there  is  no  published  paper  proposing  such 
model  and  method  for  comparison  purposes. 

However,  the  proposed  e-constraint  method  has  been  imple¬ 
mented  on  a  test  case  in  Refs.  [20,29-31  ]  comprising  4  hydro  units 
and  3  thermal  units.  The  Economic  Load  Scheduling  (ELS)  pre¬ 
sented  in  Refs.  [29-31]  only  takes  into  consideration  the  economic 
issues.  In  addition,  the  Economic  Emission  Scheduling  strategy  is 
represented  via  EES  intended  to  minimize  the  emission  issues 
while  CEES  i.e.  Combined  Economic  Emission  Scheduling  strategy 
optimizes  two  objective  functions  comprehensively  by  converting 
the  original  problem  into  a  optimization  problem  with  only  one 
objective  function  employing  a  set  of  weights.  The  nonlinear 
formulation  proposed  in  Refs.  [20,29-31]  is  used  to  make  a 
reasonable  comparison  with  the  results  obtained  in  these  works. 

The  optimal  solution  obtained  from  e-constraint  method  con¬ 
sidering  equal  weighing  factors  by  DM  is  depicted  in  Fig.  10.  The 
power  generated  by  the  first  hydro  unit  and  the  first  thermal  unit 
for  each  hour  is  denoted  by  Phlt  and  Pslt,  respectively.  The  ELS 
and  EES  can  be  calculated  using  the  proposed  method  by  optimiz¬ 
ing  economic  load  scheduling  and  economic  emission  scheduling, 
respectively.  Note  that,  the  payoff  table  is  constructed  by  the 
results  obtained  from  ELS  and  EES.  The  results  of  the  proposed 
method  are  represented  in  Table  4  showing  the  superiority  in  the 
case  of  quantity  for  all  cases  compared  to  those  obtained  in 
Refs.  [20,29-31],  For  instance,  the  costs  obtained  using  this 
approach  is  7051.72  (47906-40854.28  =  7051.72),  4059.72,  2652.72 
and  3489.72($)  which  are  less  than  the  ones  obtained  in 
Refs.  [20,29-31],  Furthermore,  the  emission  resulted  via  this  method 
is  less  than  the  ones  reported  in  Refs.  [20,29-31  ]  i.e.  10053.8,  3434.8, 
2002.8  and  1227.8  (lbs).  Also,  the  proposed  method  is  better  than 
other  methods  in  the  case  of  solution  time.  The  optimization  statistics 
for  multi-objective  SCUC  problem  on  this  test  case  are  represented  in 
Table  6.  The  results  obtained  from  implementing  the  proposed 
approach  demonstrate  that  the  optimal  daily  generation  scheduling 


of  hydrothermal  systems  can  be  efficiently  solved  using  lexicographic 
optimization  and  augmented  e-constraint  method  while  the  solution 
time  is  rationale. 


5.  Conclusion 

The  lexicographic  optimization  along  with  augmented 
e-constraint  technique  are  proposed  in  this  paper  and  prosper¬ 
ously  utilized  to  solve  the  multi-objective  SCUC  problem.  Also, 
simulations  have  been  done  on  the  IEEE  188-bus  system  with  54 
thermal  units  and  8  hydro  plants  to  verify  the  effectiveness  of  the 
proposed  method.  The  objective  of  ISO  is  to  compromise  the 
conflicting  objectives  of  cost  minimization  and  gaseous  emissions 
minimization.  The  multi-objective  SCUC  problem  is  effectively 
solved  employing  the  proposed  method  and  Pareto  optimal  solu¬ 
tions  are  generated.  Afterwards,  a  Fuzzy-based  decision  making 
procedure  is  implemented  to  select  the  most  desired  solution 
among  Pareto  optimal  solutions.  The  obtained  results  using 
lexicographic  optimization  and  augmented  E-constraint  approach, 
in  comparison  with  the  ones  obtained  from  an  interactive  fuzzy 
satisfying  method  [20],  differential  evolution  method  [29], 
quantum-behaved  particle  swarm  optimization  method  [30]  and 
cultural  algorithm  [31],  confirms  the  capability  of  the  proposed 
method  both  in  the  case  of  precision  and  execution  time.  The 
ongoing  research  work  is  to  introduce  security  indices  (overload 
index,  voltage  drop  index  and  voltage  stability  margin)  to  the 
proposed  problem. 
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